home *** CD-ROM | disk | FTP | other *** search
/ The X-Philes (2nd Revision) / The X-Philes Number 1 (1995).iso / xphiles / hp48_1 / gamma < prev    next >
Internet Message Format  |  1995-03-31  |  2KB
  1. From comp.sys.handhelds Tue Feb  5 09:10:38 1991
  2. Path: mentor.cc.purdue.edu!purdue!bu.edu!rpi!zaphod.mps.ohio-state.edu!swrinde!cs.utexas.edu!uunet!shelby!agate!dog.ee.lbl.gov!usenet
  3. From: austern@ux5.lbl.gov (Matt Austern)
  4. Newsgroups: comp.sys.handhelds
  5. Subject: Gamma function
  6. Summary: Here's a program to compute gamma(z) for arbitrary complex z.
  7. Message-ID: <9483@dog.ee.lbl.gov>
  8. Date: 4 Feb 91 06:17:59 GMT
  9. Reply-To: austern@ux5.lbl.gov (Matt Austern)
  10. Organization: Lawrence Berkeley Laboratory (theoretical physics group)
  11. Lines: 51
  12.  
  13. The gamma function is one of the standard "special functions" of
  14. complex analysis.  It is therefore a little bit frustrating to see
  15. that the gamma function on the HP 48 only accepts real arguments.
  16. (Not that I'm complaining, of course: HP is the only company I know
  17. of that puts the gamma function on its calculators at all.)
  18.  
  19. The following program corrects this.  The algorithm isn't as accurate
  20. as that on the HP 48, but it still is quite good: better than 1 part
  21. in 10^9.
  22.  
  23. This program uses an algorithm from Numerical Recipes.  (Press, et.
  24. al.) They give no theory behind the algorithm, besides the obvious
  25. observation that it asymptotically approaches Stirling's
  26. approximation, but they do give a reference.  It is only valid for 
  27. Re (z) > 1, so I use the "reflection formula" if given an argument 
  28. with Re(z) < 1.
  29.  
  30. You'll note that this program calls itself (in the reflection
  31. formula), so you should either store it under the name GAMMA or
  32. change that call.
  33.  
  34.  
  35. %%HP: T(3)A(R)F(.);
  36. \<< \-> z1
  37.   \<< z1 1 - \-> z
  38.     \<<
  39.       IF z RE 0 <
  40.       THEN \pi \->NUM z
  41. * DUP SIN / 1 z -
  42. GAMMA /
  43.       ELSE z 5.5 +
  44. z .5 + ^ z 5.5 +
  45. NEG EXP * 2 \pi \->NUM
  46. * \v/ * 1 76.18009173
  47. z 1 + / +
  48. -86.50532033 z 2 +
  49. / + 24.01409822 z 3
  50. + / + -1.231739516
  51. z 4 + / +
  52. .00120858003 z 5 +
  53. / + -.00000536382 z
  54. 6 + / + *
  55.       END
  56.     \>>
  57.   \>>
  58. \>>
  59. -- 
  60. Matthew Austern    austern@lbl.bitnet     Proverbs for paranoids, 3: If 
  61. (415) 644-2618     austern@lbl.gov       they can get you asking the wrong 
  62.                                           questions, they don't have to worry
  63.                                           about answers.
  64.  
  65.